Broome et al. BioMedical Engineering OnLine 2013, 12:69 
http://www.biomedical-engineering-online.eom/content/12/1/69 



E BioMedical Engineering 
OnLine 



RESEARCH 



Open Access 



Closed-loop real-time simulation model of 
hemodynamics and oxygen transport in the 
cardiovascular system 



Michael Broome l,2,y , Elira Maksuti" 5 , Anna Bjallmark" 5 ' 4 , Bjorn Frenckner l,!) and Birgitta Janerot-Sjoberg 



3,4 



.3,6,7 



* Correspondence: 
michael.broome@karolinska.se 
^CMO Department, Karolinska 
University Hospital, SE-171 76, 
Stockholm, Sweden 
Anesthesiology and Intensive Care, 
Department of Physiology and 
Pharmacology, Karolinska Instituted 
Stockholm, Sweden 
Full list of author information is 
available at the end of the article 



Abstract 

Background: Computer technology enables realistic simulation of cardiovascular 
physiology. The increasing number of clinical surgical and medical treatment options 
imposes a need for better understanding of patient-specific pathology and outcome 
prediction. 

Methods: A distributed lumped parameter real-time closed-loop model with 26 
vascular segments, cardiac modelling with time-varying elastance functions and 
gradually opening and closing valves, the pericardium, intrathoracic pressure, the 
atrial and ventricular septum, various pathological states and including oxygen 
transport has been developed. 

Results: Model output is pressure, volume, flow and oxygen saturation from every 
cardiac and vascular compartment. The model produces relevant clinical output and 
validation of quantitative data in normal physiology and qualitative directions in 
simulation of pathological states show good agreement with published data. 

Conclusion: The results show that it is possible to build a clinically relevant real-time 
computer simulation model of the normal adult cardiovascular system. It is 
suggested that understanding qualitative interaction between physiological 
parameters in health and disease may be improved by using the model, although 
further model development and validation is needed for quantitative patient-specific 
outcome prediction. 

Keywords: Cardiovascular simulation, Time-varying elastance functions, Lumped 
parameter model, Valve model, Oxygen transport model, Computer simulation 



Background 

The increasing number of surgical and medical treatment options in cardiovascular 
disease imposes a need for better understanding of patient-specific pathology and out- 
come prediction. Due to the complexity of the cardiovascular system, it is not always 
possible to understand the role of every cardiovascular parameter; additionally, predic- 
tion of the hemodynamic features can be even more challenging when several of those 
parameters change simultaneously. In many cases the use of simulation models can 
improve the understanding of both normal physiology and pathophysiology. Moreover, 
clinical training in simulation centres, where simulation models have been connected 
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to mannequins responding to therapeutic measures, has proved to further increase 
diagnostic and therapeutic accuracy [1]. 

Other authors [2] have published complex hemodynamic models based on a reduc- 
tionistic approach, where contractions are based on cellular properties including elec- 
trophysiology, ion channels, myofilament structure and fiber orientation. These models 
are relevant for exploring normal and disturbed contractile mechanisms as well as inte- 
grated organ function. The computational load in 3D models is however presently too 
large to allow clinical and educational use in studies of central hemodynamics. Simpler 
models based on cellular properties [3] tend to be less relevant on higher explanatory 
levels when studying integrated properties of the entire cardiovascular system, although 
being computationally less demanding. 

Models focusing on particular aspects of hemodynamics such as mitral valve dynam- 
ics [4], remodeling of heart and vessels [5], individual drug effects [6] and effects of 
therapeutic measures such as positive end-expiratory pressure (PEEP) [7] offer good ex- 
amples of how cardiovascular simulation models can be validated in specific situations 
and used clinically. It is however difficult to know to what extent these models can be 
used in clinical situations other than those described in the articles. 

The aim of this paper is to present a closed-loop, real-time, lumped parameter car- 
diovascular simulation model, which includes the behavior of the four cardiac cham- 
bers, interatrial and interventricular septum, pericardium, intrathoracic pressure 
changes, heart valves, intracardiac shunts and vascular system. The current model is 
based on previous publications [8,9]. The circulatory system interacts with pericardial 
as well as intrathoracic pressures as in one of the previous models [9] and several 
improvements have been made. The vascular system has been expanded and the circu- 
latory system closed in order to reach more realistic steady state output when parame- 
ters in the model are changed. The systolic cardiac model is replaced in order to more 
realistically simulate cardiac work and valvular timing. Passive diastolic properties are 
now non-linear [10] resulting in more realistic changes of filling pressures with ven- 
tricular dilatation. Atrial septal interactions are included in addition to a modified ven- 
tricular septum [9]. The valves are gradually opening and closing [11] and adjustment 
of valvular properties therefore allow more realistic simulation of both stenotic and 
regurgitant valves. Potential atrial septal defects, ventricular septal defects and persist- 
ent ductus arteriosus have been included. An oxygen transport model allowing studies 
of oxygen delivery and uptake including myocardial oxygen balance have been devel- 
oped. The most unique feature of this model is however not any particular advanced 
modeling solutions in each sub-part, but rather the holistic approach giving a compre- 
hensive overview of cardiovascular hemodynamics in health and disease. The possibility 
to change the model parameters one by one facilitate understanding in educational 
sessions, but real-time simulation and fast feed-back after parameter changes are also 
important features in a clinical decision support scenario in cardiac surgery and 
intensive care. 

The model is presented from a normal cardiovascular physiological point of view, but 
a limited number of pathological scenarios will also be shown. A secondary aim is to 
show that the model generates valid output, which is consistent with earlier published 
clinical and experimental data, focusing on the behavior of the entire system rather 
than on single variables. 
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Methods 

General overview of the model 

The model (Figure 1) is based on the following compartments: the four cardiac chambers 
and corresponding valves, the pericardium, the systemic circulation consisting of the aortic 
root, ascending aorta, proximal aortic arch, right carotid/subclavian artery, distal aortic arch, 
left carotid/subclavian artery, descending aorta and the representation of a peripheral arter- 
ial system (See also electrical analogue in Additional file 1: Figure SI). The peripheral arterial 
system contains resistance arteries (arteriolae), capillaries and small veins emptying into the 
inferior caval vein, while the two carotid/subclavian arteries in a similar way end up in the 
superior caval vein. The pulmonary circulation is represented by the main pulmonary trunk, 
pulmonary resistance arteries (arteriolae), capillaries, small veins and pulmonary veins 
emptying into the left atrium. Two coronary vessels are also included. Input to the vascular 
model are physical dimensions [12,13], tissue properties and number of parallel vessels in 
each compartment modified to fit the anatomical structure of the present vascular model; 
while pressure, flow, volume and oxygen saturation are simulated outputs and can be 
displayed. The cardiac, valvular, pericardial and vascular properties can be modified in detail 
by changing the parameters presented in the supplement (Additional file 2: Table Sl-5). 
Model complexity is chosen to enable future simulations of central hemodynamics in major 
congenital and acquired circulatory pathologies and extra-corporeal circulatory support. 



Cardiac model 
The heart chambers 

The four cardiac chambers passive and active properties are modeled with separate 
time-varying elastance functions based on the "Double-Hill equation" as previously 
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Figure 1 Sketch showing the cardiac and vascular components of the simulation model. The dark 
yellow area is the pericardium containing the cardiac chambers and coronary vessels. The light yellow area 
is the intra-thoracic space containing the pericardium, the pulmonary circulation and the thoracic aorta. The 
extra-thoracic space contains the carotid/subclavian circulation in the upper part and the rest of the 
systemic circulation in the lower left part. 
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described [11,14] including a slightly modified Starling mechanism (Figure 2) [8]. The 
model includes internal chamber flow resistance and viscous wall properties of the 
heart, in accordance with previous models [9,11]. Details are shown in the cardiac 
model section of the supplement. 

By using a continuous function through the cardiac cycle manual setting of valvular 
timing [9] is avoided. Each cardiac chamber has its own elastance function of realistic 
shape (Additional file 3: Figure S3). Timing and duration of contractions in the cardiac 
chambers can be set independently, in order to create variable atrioventricular (AV) 
and interventricular ( W) intervals as well as realistic timing of valve opening and closure. 
Thereby the model offers additional opportunities to simulate arrhythmias as well as 
interventricular dyssynchronicity. As an example atrial fibrillation can be simulated by 
turning off atrial contractions and randomizing the time between ventricular contractions. 

The electrocardiogram (ECG) (Additional file 4: Figure S4; Additional file 5: Figure 
S5) is simulated by summing the four cardiac elastance functions with a linear weight 
factor, based on the approximate distance from each cardiac chamber to the chosen 
surface ECG lead assuming an electromechanical delay of 60 ms. The ECG does not by 
itself influence the simulated physiology and is included as an output of the model only 
to add clinical realism during educational sessions, therefore the detailed algorithms 
are not described. 

The heart valves 

Valve pressure gradients are composed of a Bernoulli resistance and an inertial term 
(Additional file 2: Equation S4 and Figure S6 of Additional file 6) as earlier described 
[8,9]. Minimum and maximum valvular areas can be set to simulate valvular stenosis 
and regurgitation. The valve areas change gradually (Additional file 4: Figure S4) as a 
function of pressure gradients and valve inertial constants [11] as is shown in the 
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Figure 2 Left ventricular pressure-volume loops during a preload reduction maneuver in a case 
with slightly reduced left ventricular function (E max 1.4 and E min 0.08). The diastolic and end-systolic 
pressure-volume relations are shown in thin black lines. The end-systolic line is curved due to the "heart law 
of Starling" reducing contractile strength during volume overload according to Additional file 2: Equation 
S2. The curved end-diastolic relation illustrates the increase in passive stiffness when the ventricle is dilated 
according to Additional file 2: Equation S3. 
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supplement (Additional file 2: Table S4). Blood flow inertance L is determined by dens- 
ity p, the variable inflow length / and the valve area A (Equation 1). Inflow length / is 
set to a value identical to the instantaneous diameter of the valve in contrast to the 
constant value in Mynard et al [11]. 

L= P_l = P_2^ = P_2 

A TIT 1 71 r 



The atrial and ventricular septum 

The interventricular septal displacement is simulated using pressure transmission 
through the muscular septum as previously described [9,15] (Equation 2). As a further 
development the ventricular septal stiffness E sv is modeled in order to increase in pro- 
portion to the left ventricular systolic elastance e iv (t) (Equation 3), thereby stiffening 
during contraction, in accordance with known physiological characteristics. Additional 
file 5: Figure S5 shows the ventricular septal shift at different stiffness. 

Esv e lv / n \ 

p * = n^r' eb, ' viv+ n^r' p "' (2) 

J-'SV i c /v ^sv i c /v 

E sv = E SV Q-ei v (t) (3) 

Atrial septal interaction is modeled in a similar way. The septal interaction does not 
give any significant contribution during normal physiology in apnea, but is important 
when simulating interventricular dyssynchrony, variations in intrathoracic pressures in- 
cluding artificial ventilation and changes in ventricular loading conditions. 

The pericardium 

In order to resemble pericardial function in health and disease a user-defined exponen- 
tial function relating intra-pericardial pressure to total heart volume is adopted from 
the literature [9] (Additional file 7: Figure S7). As a further development of the previous 
model the minimal pericardial pressure is allowed to be a negative value in agreement 
with what is found experimentally when measured during intrathoracic pressure 
changes or hypovolemia [16]. Further details can be found in the supplement. 



Intracardiac and arteriovenous shunts 

Atrial septal defects, ventricular septal defects and persistent ductus arteriosus flows are 
calculated according to Equation 4, with instantaneous flow Q proportional to the area A 
of the defect, the square-root of the pressure difference AP and the Gorlin empirical con- 
stant C assumed to be 1 [17]. Acceleration g due to gravity is set to 980 cm/s 2 . 

Q = C-A-y/2-g-AP = 4A3A-VAP (4) 



Vascular model 

The systemic circulation 

Transmural pressure p in each vascular compartment is related to the variable 
segmental vascular volume v according to an exponential relation (Additional file 2: 
Equation S7, Additional file 8: Figure S12) [9] creating an increase in stiffness with pro- 
gressive distension of the vascular wall. 
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Resistance R 0 , inertia I 0 , and elastance E 0 at normal mean pressure P 0 are calculated 
from published and estimated vascular properties [12,13], and relations [18] described 
in Equations 5-7. These vascular properties are all approximately inversely propor- 
tional to n, which is the number of parallel vessels, lumped into each compartment 
[18]. The parameter rj is blood viscosity, /, r 0 and h is the length, radius and thickness 
of a single vessel respectively and p is blood density. 



71 ' Tq H 



Y inc h 



2 tt Tq I n 



(7) 



The same value of Youngs modulus Y inc is used in every vascular segment, implying 
that vascular stiffness properties are solely explained by this constant and vascular di- 
mensions. The actual volume-dependent segmental elastances, resistances and inertias 
are updated in each calculation step based on the volume and radius assuming constant 
vessel length and thickness. Details can be found in the blood vessels section of the 
supplement (Additional file 2: Table S5-7). 

The resistance £1, in series with the capacitor situated in parallel to the blood flow in 
every vascular segment, corresponds to the viscous property of the vascular wall, damp- 
ening the flow pulse. The parameter value assigned to this component is calculated as 
in Equation 8 with addition of a scaling factor A as compared to the characteristic im- 
pedance in a classical 3-component Windkessel model [19,20]. 

cl = \Vl 7 e (8) 

The same formula was used to calculate Cl in every vascular segment as previously sug- 
gested [19], although the name characteristic impedance is less well suited in a distributed 
model. A scaling factor X between 0.0 and 1.0 is used to tune the vascular damping prop- 
erties. The effects of changing the scaling factor X can be seen in Figure 3. 

The systemic capillaries and veins are described by the same principal vascular model 
as the arteries although inertial effects are of less importance due to lower pulsatility 
[21] and a larger cross-sectional area. Every single vessel has high resistance and stiff- 
ness; however the large amounts of vessels in parallel lumped into each compartment 
result in both a low resistance (Equation 5) and a low elastance (Equation 7), 




A=0.02 A=0.10 A=0.50 A=1.00 

Figure 3 Effects of changing vascular damping factor A on cardiac and vascular pressures 
(peripheral artery (light red), ascending aorta (dark red), left ventricle (red), pulmonary artery 
(orange), right ventricle (yellow), left atrium (brown), right atrium (blue)) related to the viscous 
properties of the vascular walls. Very low values results in numerical instability in the model. A value of 
0.5 is chosen in all presented simulations. 
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corresponding to high compliance and being related to their physiological volume res- 
ervoir function. The distribution of the blood volume between arteries, capillaries and 
veins in the systemic and pulmonary circulation in the model is shown in Additional 
file 2: Table S7. 

The pulmonary circulation 

The pulmonary artery, resistance vessels, capillaries and veins are modeled as in the sys- 
temic circulation. An adjustable pulmonary shunt mimics the circulation without gas 
exchange due to anatomical anomalies or ventilation-perfusion disturbances as in 
atelectatic lung tissue. 

The coronary circulation 

Coronary circulation is simulated with a left and right coronary artery emptying in the 
right atrium as described in the supplement (See also Additional file 9: Figure S8). 
Myocardial oxygen consumption is calculated according to the pressure-volume area 
concept of Suga et al. [22] on a beat-by-beat basis and estimations of myocardial oxy- 
gen balance can therefore be calculated based on both supply and demand. 

Oxygen transport 

Blood oxygen content in every compartment is calculated as a constant multiplied by 
the hemoglobin (Hb), the oxygen saturation (Sat) and the volume of the compartment 
(Volume) as in clinical routine and most catheterization laboratories [23] (Equation 9). 

C = 0.0000136 • Hb ■ Sat ■ Volume (9) 

Oxygen saturation within each compartment is considered homogenous and oxygen 
transport between compartments is proportional to blood flow, hemoglobin level and 
saturation. Calculations concerning blood mixing are simplified as effects of blood jets, 
incomplete mixing and physically dissolved oxygen not are taken into account. 

Baroreceptor reflex 

A baroreceptor reflex can be activated in the model that affects heart rate, cardiac con- 
tractility (maximum elastance) and arterial vascular resistance as in Sun et al. [9] and 
described in the baroreceptor section of the supplement. 

Simulation 

The simulator was compiled as a stand-alone software developed in Visual Basic.NET 
2010 (Microsoft Corporation, Redmond, WA, USA) and .NET Framework 4.0 
(Microsoft Corporation, Redmond, WA, USA). A time step of 0.25 ms was used to as- 
sure stability. The implemented 62 differential equations were solved with an implicit 
Euler numerical method [24]. Mean values were calculated as running arithmetical 
means to simplify calculations and conserve memory resources. The program version 
used was Aplysia CorVascSim 4.4.0.62 (Aplysia Medical AB, Stockholm, Sweden). The 
model can be run in real-time on a standard personal computer (example: Intel pro- 
cessor Core i3 3.06 GHz, RAM 3 Gb, Graphics ATI Radeon HD4550, Windows 7 32/ 
64-bit). 
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Parameter settings in simulations of normal physiology, systolic heart failure, diastolic 
heart failure, aortic stenosis, aortic regurgitation, a Valsalva maneuver, exercise and 
progressive arteriosclerosis 

All hemodynamic data were presented with zero intrathoracic pressure in order to resem- 
ble end-expiratory values. With the basic set of parameters presented in the supplement 
an adult man weighing 70 kg is described with a blood volume of 5600 ml (80 ml/kg), 
blood viscosity of 0.00024 mmHg-s and blood density of 1.060 g/cm 3 . Systolic left 
ventricular heart failure was simulated by decreasing the maximum elastance from 2.8 to 
1.0 mmHg/ml. Diastolic left ventricular heart failure was simulated by increasing the basal 
passive diastolic elastance E min from 0.05 to 0.12 mmHg/ml. Aortic stenosis was simu- 
lated by decreasing the open aortic valve area from 5.0 to 0.7 cm 2 . Aortic regurgitation 
was simulated by increasing the closed aortic valve area from 0.0 to 0.2 cm 2 . Simulation 
of a Valsalva maneuver was accomplished by increasing intrathoracic pressure from 0 to 
10 mmHg for approximately 20 seconds with and without baroreceptors activated. All 
other parameters and settings were identical. 

The stepwise simulation of exercise was accomplished by increasing heart rate from 
72 to 144 min" 1 , increased left and right ventricular contractility by 50% and increased 
systemic and pulmonary resistance arterial diameters by 50%. The cardiac elastance pa- 
rameters a 1 and a 2 were scaled with the square-root of the factor changing the heart 
rate to preserve a realistic relation between systolic and diastolic times. The vascular 
damping scaling factor A was set to 0.5 in all vascular segments and presented simula- 
tions. The oxygen saturation in the pulmonary capillaries was set to 99.4%, correspond- 
ing to a normal alveolar p02 of 13.3 kPa and the pulmonary shunt was set to 10% of 
cardiac output in all presented simulations. Hemoglobin level was set to 140 g/1 and 
systemic oxygen consumption to 250 ml/min in addition to myocardial oxygen extrac- 
tion. A Youngs modulus Y inc of 3000 mmHg (= 0.40 MPa) similar to the value in 
Wang et al. [12] was used to calculate the vascular elastance E 0 at the vessel-specific 
normal mean pressure P 0 (see Additional file 2: Table S5) in all simulations above. A 
progressive arteriosclerotic process was simulated by changing vascular stiffness in all 
vascular segments through variations in Youngs modulus Y inc between 2000 and 
6000 mmHg. Note that comparisons with real clinical cases will not be performed due 
to unknown patient-specific secondary adaptations. 

Sensitivity analysis 

The relation between changes in input parameters x and model output variables y was 
studied with a sensitivity analysis during simulation of normal physiology (Additional 
file 10). Input parameters were increased 10% one by one. Seven model output variables 
[left ventricular end-systolic pressure (LVESP), right ventricular end-systolic pressure 
(RVESP), left atrial pressure (LAP), right atrial pressure (RAP), left ventricular stroke 
work (LVSW), right ventricular stroke work (RVSW) and cardiac output (CO)] focusing 
on cardiac function were evaluated after 60 seconds assuring a hemodynamic steady- 
state. Sensitivity S was calculated according to Equation 10. 




(10) 
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Results 

Overall, the direction and magnitude of the changes in output variables and curve 
shapes were within a plausible clinical range in the simulated cases. The specific data 
for each case is presented below. 



Normal physiology 

Running the model with basal settings representative of a 30 years old healthy man 
resulted in cardiac output 5.09 1/min, systemic arterial blood pressure 112/61(79) mmHg, 
pulmonary arterial pressure 24/8(11) mmHg, mean right atrial pressure 4 mmHg and 
mean left atrial pressure 5 mmHg. The model output of pressures (Figure 4 and 
Additional file 11: Figure S9), volume-flow rates (Figure 5 and Additional file 12: Figure S10) 
and timing within the cardiac chambers and the systemic and pulmonary circulatory parts 
created by the model resembled in pattern and magnitude the ones found in clinical 
practise (Table 1). 

The shapes and areas (stroke work) of the ventricular pressure-volume loops 
(Additional file 13: Figure SI la) as well as right and left ventricular ejection fraction were 
within normal limits [25] (Table 1). The atrial loops (Additional file 13: Figure SI lb) were 
similar to the scarce published data [26]. 

Systolic heart failure 

The simulation of systolic heart failure resulted in decreased left ventricular ejection 
fraction, stroke work, stroke volume and systolic pressure, while left ventricular end- 
diastolic volume and filling pressures increased (Figure 6 and Table 2). 
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Figure 4 Left ventricular (red), aortic root (pink) and left atrial pressure (brown) during two heart 
cycles. An early diastolic as well as a late positive ventriculo-atrial pressure gradient results in the mitral 
E and A-waves respectively. A post-ejection aortic pressure incisura is seen corresponding to aortic valve 
closure followed by a "reflected" wave. 
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Figure 5 Aortic valve flow (red), mitral valve flow (brown) and pulmonary vein flow (purple). Small 

regurgitant flows are seen when the valves are closing. A minimal reversal of pulmonary vein flow is seen 

corresponding to the left atrial contraction (mitral A-wave). 
\ J 



Diastolic heart failure 

The simulation of diastolic heart failure showed preserved left ventricular ejection frac- 
tion, but decreased stroke work, end-diastolic volume, stroke volume and systolic pres- 
sure, while left ventricular filling pressures increased (Figure 7 and Table 2). 



Aortic stenosis 

The simulation of aortic stenosis resulted in decreased left ventricular ejection fraction, 
stroke volume and systemic blood pressure, although left intraventricular pressure and 
stroke work increased (pressure-load). Left ventricular end-diastolic volume and filling 
pressures increased (Figure 8 and Table 2). 



Aortic regurgitation 

Simulated data of aortic regurgitation showed apparent left ventricular ejection fraction 
and stroke volume increase, systemic blood pressure decrease although stroke work in- 
creased (volume-load). Left ventricular end-diastolic volume and filling pressures in- 
creased (Figure 8 and Table 2). 



Valsalva 

The simulation of a Valsalva maneuver without and with activated baroreceptor reflex 
is shown in Figure 9. The sudden decrease in afterload with preserved preload when in- 
creasing intrathoracic pressure [27] is the explanation for the initial rise in arterial pres- 
sure as indicated in the figure. The opposite effect was seen immediately after decreasing 
intrathoracic pressure. Arterial blood pressure was as expected better preserved when the 
baroreceptor reflex was activated due to increase in heart rate, ventricular contractility 
and systemic arterial resistance. 
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Table 1 Results from a simulation of a healthy man (age 30 years, weight 70 kg, 170 cm) 
compared to published normal values [13,25] 


Global 


Unit 


Normal values 


Simulated normal case 


Heart rate 


min" 1 


69 ±17 


72 


Cardiac output 


I • min" 1 


6.3 ± 2.4 


5.09 


Stroke volume 


ml 


89 ±30 


71.6 


Systolic arterial pressure 


mmHg 


112/74 (87) 


112/61 (79) 


Pulmonary arterial pressure 


mmHg 


22/9 (14) 


24/8 (11) 


Myocardial volume 


ml 


374+ 1 10 


380 


Maximum total heart volume 


ml 


762 ±133 


737 


Total blood volume 


ml 


5105 ± 578 


5600 


Systemic arterial volume 


% 


11.4 


7.0 


Systemic capillary volume 


% 


5.4 


4.4 


Systemic venous volume 


% 


70.0 


71.9 


Pulmonary blood volume 


% 


8.1 


10.7 


Pulmonary capillary volume 


% 


2.5 


2.2 


Systemic arterial 0 2 saturation 


% 


96.4 ±0.5 


96.2 


Mixed venous 0 2 saturation 


% 


70 ±5 


68.3 


Coronary sinus 0 2 saturation 


% 


32 ±8 


28.6 


Left ventricle 








Ejection fraction 


% 


59 ±4 


67 


Stroke work 


mmHg-ml 


5600 ±1000 


6735 


End-systolic elastance 


mmHg-ml" 1 


3.51 ±1.26 


2.43 


Preload recruitable stroke work 


mmHg 


123 ±36 


88.8 


dp/dt max 


mmHg-s" 1 


1840 ±327 


2359 


dp/dt min 


mmHg-s" 1 


-1864 ±390 


-2666 


Tau 


ms 


33 ±8 


18 


End-diastolic volume 


ml 


123 ±28 


106.8 


LA volume min — max 


ml 


29.9 - 76.7 


55 - 81 


LAP 


mmHg 


8.3 ± 2.5 


4.8 


LVTei 


— 


0.39 ± 0.05 


0.41 


LV Ees/Ea 


— 


1 .62 ± 0.80 


1.78 


Right ventricle 








RV EF 


% 


61.0 ±5.8 


68 


RVSW 


mmHg-ml 


1980 


1706 


RV E es 


mmHg-ml" 1 


0.66 


0.51 


RV dpdt max 


mmHg-s" 1 


422 ± 46 


538 


RV dpdt min 


mmHg-s" 1 


-355 ±45 


-289 


RV Tau 


ms 




33 


RV EDV 


ml 


111 ±22 


105.6 


RA volume min — max 


ml 




49 - 78 


RAP 


mmHg 


4.8 ± 2.4 


3.7 


RVTei 




0.28 ± 0.05 


0.13 


RV Ees/Ea 






1.14 



Values are presented as mean ± standard deviation when available. 
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Figure 6 Intracardiac pressure-volume loops from the left ventricle (left panel) and left atrium 
(right panel) in simulation of pure systolic left heart failure (E max 2.8— ► 1.0). Stroke volume and blood 
pressure decrease. Left ventricular and atrial volumes as well as filling pressures increase. 



Exercise 

The stepwise simulation of exercise is shown in Figure 10. Starting from a resting state, 
heart rate and cardiac contractility is increased followed by reductions in systemic and 
pulmonary vascular resistance resulting in a total cardiac output above 10 1/min. 



Arteriosclerosis 

The presence of arteriosclerosis, simulated by increased vascular stiffness, resulted in 
increased systolic aortic pressure and decreased diastolic aortic pressure. A decrease in 
systolic ejection time during this process is indicated by the premature arrival of the 
reflected pressure wave changing the shape of the pressure curves in a characteristic 
way (Figure 11). 

Sensitivity analysis 

Blood volume was the single most important parameter in the sensitivity analysis 
(Additional file 10). The pericardial volume constant v_pcO was second most important. 
Pulmonary arteriolae radius (radius_pur) was more important than systemic arteriolae 
radius (radius_r). Right ventricular passive diastolic stiffness (eO_rv) and systolic con- 
tractility (emax_rv) was more important than their left sided equivalents (eO_lv and 
emax_lv, respectively). Parameters with less than 10% sensitivity are not shown. 

Discussion 

This work presents a closed-loop real-time lumped parameter simulation model of the 
cardiovascular system where normal physiology and oxygen transport as well as systolic 
heart failure, diastolic heart failure, valve failure, intracardiac shunts and pulmonary 
failure can be simulated separately or in combination. 

The validation and appreciation of clinical relevance of the model is based on empir- 
ically relevant changes in simulated clinical parameters when being compared with 
measurements made by clinically available methods such as intravascular pressure cath- 
eters and ultrasound. It is not possible to validate the large number of parameters 
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Table 2 Key observations during simulation of pathological states 



Case 


HR 


SAP 


LAP 


LVEF 


LVSV 


LVSW 


LVEDV 


LVESP 




min" 1 


mmHg 


mmHg 


% 


ml 


mmHg-ml 


ml 


mmHg 


Normal physiology 


72 


112/61(79) 


4.8 


67 


71 


6735 


107 


94 


Systolic heart failure 


72 


84/51(62) 


8.1 


39 


52 


3499 


133 


77 


Diastolic heart failure 


72 


90/53(66) 


9.2 


69 


66 


4194 


83 


72 


Aortic stenosis 


72 


89/55(68) 


7.1 


46 


58 


9117 


126 


181 


Aortic regurgitation 


72 


110/42(65) 


7.1 


76 


101 


8842 


133 


83 



HR heart rate, SAP systemic arterial pressure, LAP left atrial pressure, LVSV left ventricular stroke volume, LVSW left 
ventricular stroke work, LVEDV left ventricular end diastolic volume, LVESP left ventricular end systolic pressure. 



unambiguously with available clinical and/or experimental data, but when taking the 
complexity of the system into account and the variety of disease states that can be sim- 
ulated, the clinical relevance becomes a powerful validation and hypothesis creating 
tool. Another approach to challenge the adequacy of a cardiovascular simulation model 
is to study if model responses are adequate when evaluating model output in new sce- 
narios with pathological states not previously tested during construction and tuning of 
the model. We have chosen to present model output illuminating a range of normal 
and pathophysiological states, illustrating the wide scope of the model. The normal 
resting state, the Valsalva maneuver and physiological changes during exercise are parts 
of normal physiology that should be handled correctly by a multi-purpose cardiovascu- 
lar simulation model. The pathophysiological states chosen represent important clinical 
problems with various backgrounds. Notable is that very few parameter changes are 
needed to simulate the basic features of these states. Further parameter changes are 
needed to simulate secondary compensatory adaptations such as e.g. ventricular hyper- 
trophy in aortic stenosis. This can be automatized in a simulation model [5], but 
manual stepwise adaptation often enables a better understanding of the relation be- 
tween the primary problem and its consequences. The lack of secondary changes such 
as blood volume increase and systemic vasoconstriction usually seen in clinical cases 
with low cardiac output can also explain the small increase in filling pressures seen in 
simulation of heart failure and aortic valve failure in the present work. 
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Figure 7 Intracardiac pressure-volume loops from the left ventricle (left panel) and left atrium 
(right panel) in simulation of diastolic heart failure (E min 0.05 — > 0.12). Stroke volume and blood 
pressure decrease. Left ventricular volume decreases, while left atrial volume and filling pressures increase. 
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Figure 8 Intracardiac pressure-volume loops from the left ventricle in simulation of aortic stenosis 
(open area 5.0 — > 0.7 cm 2 , left panel) and aortic regurgitation (closed area 0.0 — > 0.2 cm 2 , right 
panel. Stroke volume decreases in aortic stenosis despite high intraventricular pressure. Apparent stroke 
volume increases and systemic blood pressure decreases in aortic regurgitation. Left ventricular end- 
diastolic pressure and volume increase in both scenarios. 



The combined effect of elastic recoil in large elastic arteries and reflection of pressure 
waves in the peripheral vasculature is of importance in coronary physiology as well as 
in determining left ventricular afterload. A distributed lumped parameter model as 
presented in this article is usually not considered adequate for describing wave trans- 
mission and reflection in the vascular tree, although the use of an inductor/inertia in 
each vascular compartment and a high resistance in peripheral components create 
pressure and flow curves with realistic amplitudes and time delays resembling effects of 
reflected waves (Figures 4, 5 and 11). The advantage with our approach is less demand- 
ing calculations enabling real-time simulation and also the possibility to connect artifi- 
cial circulatory support to any node in the model as will be explored in future works. 

Aortic characteristic impedance is originally defined as a single purely resistive com- 
ponent describing the properties of the entire vascular tree in the high frequency do- 
main, when seen from the aortic valve [20]. Placing this resistor in line with the blood 
flow, as is usually done in the 3-compartment Windkessel model, gives adequate results 
at high frequencies, but obscures the meaning of the more basic concept "total periph- 
eral resistance" in the low frequency domain [19]. When placing the resistor in line 
with the capacitor in parallel to the blood flow in each vascular compartment as is done 
in the present model (Additional file 14: Figure S2), this component still fulfills the role 
of the characteristic impedance in the high-frequency domain without affecting the 
total peripheral resistance. This is the reason why the characteristic impedance is calcu- 
lated according to Equation 8. Furthermore it has a physical meaning representing the 
damping effect imposed by the viscous properties of the vascular wall (non-linear Max- 
well viscoelastance model), although the detailed relation to tissue properties remains 
to be elucidated. 

The algorithms for scaling of the cardiac elastance functions with changing heart rate 
seem efficient up to a frequency of 160 min" 1 . Above this heart rate cardiac output 
decreases suggesting either suboptimal diastolic heart function (in the model) or sub- 
optimal scaling of the elastance functions. It is also possible that atrioventricular (AV)- 
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Figure 9 Simulation of Valsalva maneuver (increase in intrathoracic pressure from zero to +10 
mmHg between arrows) without (left panel) and with (right panel) baroreceptor reflex activated. 

ECG (blue), arterial pressure (red), pulmonary arterial pressure (orange) and central venous pressure (blue), 
pressure increase when increasing intrathoracic pressure (yellow arrow), decrease when releasing pressure 
(red arrow). The baroreceptor reflex increase heart rate, left ventricular contractility and systemic arterial 
resistance. Blood pressure is better maintained and a pressure overshoot is seen after release of 
the pressure. 

V J 



plane motion, not included in the present model, is necessary to improve ventricular 
filling during tachycardia in real life physiology. 

Sensitivity analysis in a non-linear model should ideally be performed for every pos- 
sible set of input parameters. With many input parameters it is neither feasible to per- 
form the calculations nor to interpret them in every possible model state. The choice 
of relevant parameters therefore must rely on physical and physiological model rela- 
tions; however the relevance of the model is supported by the performed sensitivity 
analysis in normal physiology, where parameters determining cardiac filling, afterload 
and ventricular systolic performance are most important for cardiac flow and pressure 
generation (Additional file 10). 



Limitations and future development of the current model 

The scientific, clinical and educational value of the model is the possibility to study 
combinations of physiological and pathophysiological cardiac and vascular properties as 
well as the effects of potential modifications or treatments. This holistic approach en- 
ables hypothesis generation regarding the direction of changes and qualitative relations 
between parameters in a wide variety of normal and pathological states suggesting pos- 
sible usefulness as a clinical decision support tool. Detailed quantitative conclusions are 
however not justified without further validating studies. It must also be kept in mind 
that differences in the magnitude of changes between model output and patients will 
persist due to model simplifications and non-modelled parameters. Patient related pre- 
dictions or hypotheses created by any simulation model have to be weighed against 
other sources of knowledge such as clinical studies and experience. 

The mathematical functions creating the time-varying elastance functions are not 
based directly on cellular or mechanical properties of the myocardium, but should ra- 
ther be seen as a convenient and well-established way to mimic a realistic physiological 
ventricular filling and contraction pattern during variable loading conditions in a com- 
puter model. The chosen functions are scalable with changes in cardiac mechanical 
properties and heart rate. 
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(See figure on previous page.) 

Figure 10 a-e. Pressure volume loops in left (black) and right (gray) ventricle during stepwise 
simulation of moderate exercise, (a). Normal resting state. HR 72/min. CO 5.1 1 l/mint.SAP1 12/61(79). PAP 
24/8(12). (b). Increase in heart rate +100%.HR 144/min. CO 7.08 l/mint.SAP1 32/86(1 02). PAP 26/9(14). (c). 
Increase in right and left ventricular systolic function +50%. HR 144/min. CO 8.41 l/mint SAP1 53/95(1 1 7). 
PAP 29/9(15). (d). Decrease in systemic arterial resistance. Increase in systemic arterioli diameter +50%. HR 
144/min. CO 9.23 l/mint.SAP 86/48(60). PAP 30/9(15). (e). Decrease in pulmonary arterial resistance.lncrease 
in pulmonic arterioli diameter +50%. HR 144/min. CO 10.18 l/mint.SAP 94/50(65). PAP 27/5(13). 



Future improvements of the cardiac model aim at adjustment to known physiological 
relations, including the AV-plane movement contributions to the cardiac pumping 
mechanism [28] and scaling of the model to neonatal and pediatric size. A specific aim 
is also to develop and validate simulated hemodynamic responses in patients treated 
with extra-corporeal circulation during intensive care, such as left ventricular unloading 
during support with mechanical left ventricular assist devices and systemic oxygen de- 
livery and pulmonary hypertension in congenital cardiopulmonary anomalies. Further 
development of the model will also aim at implementing a semi-automatic adaptation 
to available patient-specific clinical data and from this enable prediction of the 
hemodynamic results from possible surgical and pharmacological interventions. 

During extreme conditions such as cardiac arrest, extreme hypovolemia or volume 
overload, when pressures are outside a normal range an improved model including 
simulation of vascular collapse and overdistension needs to be developed. Moreover, 
due to the lack of true 3D properties some pathology, e.g. myocardial ischemia and in- 
farctions with regional cardiac dysfunction, cannot be represented in detail. The same 
problem is inherent in many cardiac arrhythmias involving aberrant activation of the 
myocardium. 




0 333 667 1000 1333 1667 

Time (ms) 

Figure 11 Blood pressure in ascending aorta with changing vascular stiffness. When increasing 

(lighter) Young's elastic modulus Y systolic pressure increases and diastolic pressure decreases. The 

increased stiffness also results in a shorter systolic ejection due to a premature arterial recoil/reflection. The 

increase in afterload also results in decreasing left ventricular stroke volume despite increasing stroke work. 
\ J 
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Conclusions 

The results of the present publication show that it is possible to build a real-time com- 
puter simulation model illustrating the cardiovascular system, including oxygen trans- 
port, as evidenced by the base-line values almost entirely within reference limits of 
normal physiology and by relevant direction of changes in clinical parameters in a 
variety of pathophysiological states. The understanding of complex hemodynamic states 
may be further improved when known model parameter changes are studied and sim- 
ultaneously compared in detail by the user with clinical data. 

Additional files 



Additional file 1: Figure SI. Electrical analogue sketch showing both the heart and the vascular system. The 
dark yellow area is the pericardium containing the cardiac chambers and coronary vessels. The light yellow area is 
the intrathoracic space containing the pericardium, the pulmonary circulation and the thoracic aorta. 

Additional file 2: Closed-loop real-time simulation model of hemodynamics and oxygen transport in the 
cardiovascular system. Supplement. 

Additional file 3: Figure S3. Left (black thick) and right (gray thick) ventricular time-varying elastance functions 
during two heart cycles. Atrial elastance functions are shown with thin lines. The amplitude is closely related to 
contractile function while the volume-dependent basal level describes passive stiffness. 

Additional file 4: Figure S4. The ECG and valvular areas simulated during two normal heart beats. ECG (light 
blue) and area-changes during the heart cycle for aortic (red), pulmonary (yellow), mitral (brown), tricuspid (blue) 
valves are shown. 

Additional file 5: Figure S5. Modeled volume changes (black) due to pressure dependent ventricular septal 
volume shift from left to the right ventricle during two heart cycles. Septal stiffness values of 40, 30, 20 and 10 
mmHg/ml are shown. Lower septal stiffness increases septal shift. ECG (blue) is shown as reference. 

Additional file 6: Figure S6. Simulated pressure gradients over the mitral valve during diastole. The total mitral 
pressure gradient (thick black) is decomposed into a dominating Bernoulli gradient (grey) and an inertial 
component (dotted black). 

Additional file 7: Figure S7. Pericardial pressure-volume relations. The black curve illustrates the relation used in 
the present model with the range used by a single heart beat in orange. Pericardial pressure is slightly negative. 

Additional file 8: Figure SI 2. Pressure-volume relations in arterial vascular segments illustrating non-linear 
stiffness in the actual working pressure range of each vascular segment in accordance with Equation 10. 

Additional file 9: Figure S8. Blood flow in left (black) and right (gray) coronary artery during two heart cycles. 
Left coronary artery blood flow decreases during systole due to vascular compression caused by high left 
ventricular wall stress. 

Additional file 10: Sensitivity analysis. 

Additional file 11: Figure S9. Pressure changes simulated during two heart cycles in the aortic root (black), 
ascending aorta, proximal and distal aortic arch, descending aorta and a peripheral artery (light gray). The 
segments between the aortic root and the peripheral artery are shown in different shades of gray. 

Additional file 12: Figure S10. Flow pulsatility in pulmonary capillaries (pink) is larger than in the systemic 
capillaries (green). Pulmonary valve flow (orange) and tricuspid valve flow (blue) resembling left-sided flows are 
also shown. 

Additional file 13: Figure S11. a. Intracardiac pressure-volume loops from the left ventricle (black) and right 
ventricle (gray), b. Intracardiac pressure-volume loops from the left atrium (black) and right atrium (gray). 

Additional file 14: Figure S2. Electrical analogue sketches showing details of the heart chambers, valves and 
vascular compartments. The arrows show directions of blood flow. e(t); time-varying elastance representing systolic 
and diastolic chamber properties, R waN ; a resistance representing viscous chamber wall properties, R ou tfiow; a linear 
chamber outflow resistance, B; a non-linear Bernoulli valve resistance, l va | ve ; an inductance representing flow inertia,, 
R; a linear resistance, I; an inductance representing flow inertia, C; a non-linear capacitor representing vascular 
elasticity, Q; a non-linear resistance term representing viscous vascular wall properties. 



Abbreviations 

L: Blood flow inertance; p: Density; I: Valve inflow length; A: Area; E sv : Ventricular septal stiffness; E sv0 : Ventricular septal 
stiffness constant; ei v : Left ventricular systolic elastance; pi v : Left ventricular pressure; V| V : Left ventricular volume; 
P rv : Right ventricular pressure; Q: Flow; AP: Pressure gradient; C: Gorlin empirical constant; R 0 : Resistance at normal 
mean pressure; l 0 : Inertia at normal mean pressure; E 0 : Elastance at normal mean pressure; P 0 : Normal mean pressure; 
n: Number of parallel vessels in vascular compartment; r\: Blood viscosity; I: Vessel length; r 0 : Vessel radius at normal 
mean pressure; h: Vessel thickness; Y inc : Young's modulus; Q: Vascular wall resistance; A: Vascular wall resistance scaling 
factor; E: Elastance; Hb: Hemoglobin; Sat: Oxygen saturation; Volume: Model compartment volume; E min : Basal cardiac 
chamber passive stiffness; LVESP: Left ventricular end-systolic pressure; RVESP: Right ventricular end-systolic pressure; 
LAP: Left atrial pressure; RAP: Right atrial pressure; LVSW: Left ventricular stroke work; RVSW: Right ventricular stroke 
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work; CO: Cardiac output; S: Sensitivity; y: Output variable; x: Input parameter; BloodVolume: Blood volume; 
v_pcO: Pericardial volume constant; Radius_pur: Pulmonary arteriolae radius; Radius_r: Systemic arteriolae radius; 
eO_rv: Right ventricular passive diastolic stiffness; emax_rv: Right ventricular systolic contractility; HeartRate: Heart rate; 
Alpha2_rv: Duration of right ventricular relaxation; Radius_a: Systemic arterial radius; eOJv: Left ventricular passive 
diastolic stiffness; MeanPressure_v: Normal mean pressure in small systemic veins; Radius_puc: Pulmonary capillary 
radius; emaxjv: Left ventricular systolic contractility; Radius_rcardx: Right carotid arterioli radius; Radius_rcarsin: Left 
carotid arterioli radius; Radius_pusv: Pulmonary small vein radius; Number_pur: Pulmonary arterioli number; 
Length_pur: Pulmonary arterioli length; Length_r: Systemic arterioli length; lambda_rv: Right ventricular stiffness 
constant; Number_r: Systemic arterioli number; AV: Atrio-ventricular; dp/dt max : Maximum pressure derivative; dp/ 
dt min : Minimum pressure derivative; Tau: Ventricular relaxation constant Tau; EDV: Left ventricular end-diastolic volume; 
LA volume: Left atrial volume; Tei: Tei index = Myocardial performance index; Ees/Ea: Ventricular-vascular coupling; RA 
volume min — max: Right atrial volume; HR: Heart rate; SAP: Systolic arterial pressure; LVEF: Left ventricular ejection 
fraction; LVSV: Left ventricular stroke volume; LVEDV: Left ventricular end-diastolic volume; ECG: Electrocardiogram; 
PAP: Pulmonary arterial pressure. 
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